1:10  # a sequence from 1 to 10
##  [1]  1  2  3  4  5  6  7  8  9 10
10:1  # in the reverse order
##  [1] 10  9  8  7  6  5  4  3  2  1
$title = "About Front Matter";
$example = array(
  'language' => "php",
);


PAQUETE OneR para el aprendizaje automático

# install.packages("OneR")
 
library(OneR)
data <- optbin(iris)
model <- OneR(data, verbose = TRUE)

    Attribute    Accuracy
1 * Petal.Width  96%     
2   Petal.Length 95.33%  
3   Sepal.Length 74.67%  
4   Sepal.Width  55.33%  
---
Chosen attribute due to accuracy
and ties method (if applicable): '*'
summary(model)

Call:
OneR.data.frame(x = data, verbose = TRUE)

Rules:
If Petal.Width = (0.0976,0.791] then Species = setosa
If Petal.Width = (0.791,1.63]   then Species = versicolor
If Petal.Width = (1.63,2.5]     then Species = virginica

Accuracy:
144 of 150 instances classified correctly (96%)

Contingency table:
            Petal.Width
Species      (0.0976,0.791] (0.791,1.63] (1.63,2.5] Sum
  setosa               * 50            0          0  50
  versicolor              0         * 48          2  50
  virginica               0            4       * 46  50
  Sum                    50           52         48 150
---
Maximum in each column: '*'

Pearson's Chi-squared test:
X-squared = 266.35, df = 4, p-value < 2.2e-16
plot(model)

prediction <- predict(model, data)
eval_model(prediction, data)

Confusion matrix (absolute):
            Actual
Prediction   setosa versicolor virginica Sum
  setosa         50          0         0  50
  versicolor      0         48         4  52
  virginica       0          2        46  48
  Sum            50         50        50 150

Confusion matrix (relative):
            Actual
Prediction   setosa versicolor virginica  Sum
  setosa       0.33       0.00      0.00 0.33
  versicolor   0.00       0.32      0.03 0.35
  virginica    0.00       0.01      0.31 0.32
  Sum          0.33       0.33      0.33 1.00

Accuracy:
0.96 (144/150)

Error rate:
0.04 (6/150)

Error rate reduction (vs. base rate):
0.94 (p-value < 2.2e-16)

“Petal.Width” se identifica como el atributo con el valor predictivo más alto. Los puntos de corte de los intervalos se encuentran automáticamente (a través de la función incluida). Los resultados son tres reglas muy simples, pero precisas, para predecir las respectivas especies.optbin

Generador de contraseñas

passwords <- function(nl = 8, npw = 1, help = FALSE) {
  if (help) return("gives npw passwords with nl characters each")
  if (nl < 4) nl <- 4
  spch <- c("!", "\"", "#", "$", "%", "&", "'", "(", ")", "*", "+", ",", "-", ".", "/", ":", ";", "<", "=", ">", "?", "@", "[", "]", "^", "_", "{", "|", "}", "~")
  for(i in 1:npw) {
    pw <- c(sample(letters, 1), sample(LETTERS, 1), sample(0:9, 1), sample(spch, 1))
    pw <- c(pw, sample(c(letters, LETTERS, 0:9, spch), nl-4, replace = TRUE))
    cat(sample(pw), "\n", sep = "")
  }
}

set.seed(12)
passwords(help = TRUE)
[1] "gives npw passwords with nl characters each"
## [1] "gives npw passwords with nl characters each"
  
passwords(8)
Z6~T'b{'
passwords(14, 5)
D4Z&Q:K6x3jmT*
}25[e'rAW?@5e/
lv9lPeGK%N-<"h
r84r2$h8/l2_Mg
8qlsokbEV8h7#}

PAQUETE Ryacas para matematicas

# install.packages("Ryacas")

library(Ryacas)
# simplificación de expresiones
as_r(yac_str("Simplify(a*b*a^2/b-a^3)"))
[1] 0
# resolucion de ecuaciones
as_r(yac_str("Solve(a+x*y==z,x)"))
[1] "x==-(a-z)/y"
# expansion de expresiones
as_r(yac_str("Expand((x-2)^20)"))
expression(x^20 - 40 * x^19 + 760 * x^18 - 9120 * x^17 + 77520 * 
    x^16 - 496128 * x^15 + 2480640 * x^14 - 9922560 * x^13 + 
    32248320 * x^12 - 85995520 * x^11 + 189190144 * x^10 - 343982080 * 
    x^9 + 515973120 * x^8 - 635043840 * x^7 + 635043840 * x^6 - 
    508035072 * x^5 + 317521920 * x^4 - 149422080 * x^3 + 49807360 * 
    x^2 - 10485760 * x + 1048576)
# derivadas D
as_function <- function(expr) {
  as.function(alist(x =, eval(parse(text = expr))))
}
 
# redefine D function
D <- function(eq, order = 1) {
  yac_str(paste("D(x,", order, ")", eq))
}

Ahora, definimos la función (en este caso un polinomio 2x^3 - 3x^2 + 4x - 5simple), determinamos simbólicamente la primera y segunda derivada y trazamos todo:

xmin <- -5
xmax <- 5
 
eq <- "2*x^3 - 3*x^2 + 4*x - 5"
eq_f <- as_function(eq)
curve(eq_f, xmin, xmax, ylab = "y(x)")
abline(h = 0, lty = 2)
abline(v = 0, lty = 2)
 
D_eq <- D(eq)
D_eq
[1] "6*x^2-6*x+4"
## [1] "6*x^2-6*x+4"
 
D_eq_f <- as_function(D_eq)
curve(D_eq_f, xmin, xmax, add = TRUE, col = "blue")
 
D2_eq <- D(eq, 2)
D2_eq
[1] "12*x-6"
## [1] "12*x-6"
 
D2_eq_f <- as_function(D2_eq)
curve(D2_eq_f, xmin, xmax, add = TRUE, col = "green")

limites

# determine limits
yac_str("Limit(x,0) 1/x")
[1] "Undefined"
yac_str("Limit(x,0,Left) 1/x")
[1] "-Infinity"
yac_str("Limit(x,0,Right) 1/x")
[1] "Infinity"
# integration
yac_str("Integrate(x) Cos(x)")
[1] "Sin(x)"
yac_str("Integrate(x,a,b) Cos(x)")
[1] "Sin(b)-Sin(a)"

Como ejemplo, podemos probar en poco tiempo que la famosa aproximación es en realidad demasiado grande

yac_str("Integrate(x,0,1) x^4*(1-x)^4/(1+x^2)")
[1] "22/7-Pi"

ecuaciones diferenciasles

as_r(yac_str("OdeSolve(y' == y)"))
expression(C241 * exp(x))
as_r(yac_str("OdeSolve(y'' - 4*y == 0)"))
expression(C280 * exp(2 * x) + C284 * exp(-2 * x))

:::

And this block will be put on the right:

plot(iris[, -5])

::::




Uso de un filtro Pandoc Lua, haz click para expandir 
:::: {style="display: flex;"}

::: {}
Here is the **first** Div.

# ```{r}
# str(iris)
# ```

:::

::: {}
And this block will be put on the right:

# ```{r}
# plot(iris[, -5])
# ```

:::

::::



otro modelo https://pandoc.org/MANUAL.html#divs-and-spans

#EJECUTA CON INDEPENDENCIA DE LA VARIABLE, CUIDADO NO GUARDA LA VARIABLE X

(x <- 1+1)
## [1] 2

CODIGO VERBATIN

https://bookdown.org/yihui/rmarkdown-cookbook/multi-column.html

if (TRUE) {
  x <- 1:10
  x + 1
}
##  [1]  2  3  4  5  6  7  8  9 10 11



They went in single file, running like hounds on a strong scent, and an eager light was in their eyes. Nearly due west the broad swath of the marching Orcs tramped its ugly slot; the sweet grass of Rohan had been bruised and blackened as they passed.



Si tu navegador soporta este atributo, podrás editar este párrafo.



Upcoming Topics

For the new year, we have a great line up of articles!



Uso de un filtro Pandoc Lua, haz click para expandir 
{cat, engine.opts = list(file = "color-text.lua")}
Span = function(el)
  color = el.attributes['color']
  -- if no color attribute, return unchange
  if color == nil then return el end
    -- transform to <span style="color: red;"></span>
  if FORMAT:match 'html' then
    -- remove color attributes
    el.attributes['color'] = nil
    -- use style attribute instead
    el.attributes['style'] = 'color: ' .. color .. ';'
    -- return full span element
    return el
  elseif FORMAT:match 'latex' then
    -- remove color attributes
    el.attributes['color'] = nil
    -- encapsulate in latex code
    table.insert(
      el.content, 1,
      pandoc.RawInline('latex', '\\textcolor{'..color..'}{')
    )
    table.insert(
      el.content,
      pandoc.RawInline('latex', '}')
    )
    -- returns only span content
    return el.content
  else
    -- for other format return unchanged
    return el
  end
end

> Roses are [red and **bold**]{color="red"} and

> violets are [blue]{color="blue"}.


Roses are red and bold and

violets are blue.

I could use some help…, haz click para expandir 

public class Order
{
    public int OrderId { get; set; }
    public int CustomerId { get; set; }

    public List<int> Products { get; set; }
}
\```


Let’s see some code!, haz click para expandir 
print('Hello World!')
## Hello World!
Of course, it has to be Hello World, right?


INSERTAR FOTOS CON HTML, haz click para expandir





Hiring Manager
School of Ninja, Hacker’s University
404 Not Found Road, Undefined 0x1234, NA



hoy <- Sys.Date()

fecha<-format(hoy, format="%A %d %B %Y")

cat("La fecha del expediente es:" ,fecha)
La fecha del expediente es: domingo 20 febrero 2022



— Sara Teasdale

# install.packages("knitr")
# library(knit)
library(htmltools)
library(htmlwidgets)
library(knitr)
library(tufte)

https://albert-rapp.de/post/2021-10-16-exploratory-intro-plotly/

library(tidyverse)
library(plotly)
p <- mpg %>% 
  ggplot(aes(hwy, cty, fill = class)) +
  geom_jitter(shape = 21, size = 2, alpha = 0.5)

plotly_p <- ggplotly(p) 
plotly_p 
p_layout <- p %>% 
  ggplotly() %>% 
  layout(legend = list(
    x = 0.1, 
    y = -0.2, 
    orientation = "h"
  )) 
p_layout
set.seed(123)
jitter_hwy <- 2
jitter_cty <- 1
jittered_mpg <- mpg %>% 
  mutate(
    hwy = hwy + runif(length(hwy), -jitter_hwy, jitter_hwy),
    cty = cty + runif(length(cty), -jitter_cty, jitter_cty)
  )
plt <- jittered_mpg %>% 
  plot_ly() %>% 
  add_markers(x = ~hwy, y = ~cty, color = ~class)
plt 
dummy_dat %>% 
  mutate(percent_labels = scales::percent(percent)) %>% 
  ggplot(aes(x = group, y = percent, fill = category)) +
  geom_col() +
  geom_text(
    aes(label = percent_labels), 
    position = position_stack(vjust = 0.5),
    col = "white",
    fontface = "bold"
  ) + 
  scale_y_continuous(labels = scales::percent_format()) +
  scale_fill_brewer(palette = "Set1")
Error in mutate(., percent_labels = scales::percent(percent)): objeto 'dummy_dat' no encontrado

https://github.com/danielredondo/30diasdegraficos/blob/master/scripts/22_texto.R

# Modificación de http://www.sthda.com/upload/rquery_wordcloud.r

rquery.wordcloud <- function(x, type = c("text", "url", "file"),
                             lang = "english", excludeWords = NULL,
                             textStemming = FALSE, colorPalette = "Dark2",
                             min.freq = 3, max.words = 200) {
  library("tm")
  library("SnowballC")
  library("wordcloud")
  library("RColorBrewer")
  
    if (type[1] == "file") {
    text <- readLines(x, encoding = "UTF-8")
    text <- gsub("¿", "", text)
    text <- gsub("¡", "", text)
  }
  else if (type[1] == "url") {
    text <- html_to_text(x)
  } else if (type[1] == "text") text <- x

  # Load the text as a corpus
  docs <- Corpus(VectorSource(text))
  # Convert the text to lower case
  docs <- tm_map(docs, content_transformer(tolower))
  # Remove numbers
  docs <- tm_map(docs, removeNumbers)
  # Remove stopwords for the language
  docs <- tm_map(docs, removeWords, stopwords(lang))
  # Remove punctuations
  docs <- tm_map(docs, removePunctuation)
  # Eliminate extra white spaces
  docs <- tm_map(docs, stripWhitespace)
  # Remove your own stopwords
  if (!is.null(excludeWords)) {
    docs <- tm_map(docs, removeWords, excludeWords)
  }
  # Text stemming
  if (textStemming) docs <- tm_map(docs, stemDocument)
  # Create term-document matrix
  tdm <- TermDocumentMatrix(docs)
  m <- as.matrix(tdm)
  v <- sort(rowSums(m), decreasing = TRUE)
  d <- data.frame(word = names(v), freq = v)
  # check the color palette name
  if (!colorPalette %in% rownames(brewer.pal.info)) {
    colors <- colorPalette
  } else {
    colors <- brewer.pal(8, colorPalette)
  }
  # Plot the word cloud
  set.seed(1234)
  wordcloud(d$word, d$freq,
    min.freq = min.freq, max.words = max.words,
    random.order = FALSE, rot.per = 0.35,
    use.r.layout = FALSE, colors = colors
  )

  invisible(list(tdm = tdm, freqTable = d))
}
#++++++++++++++++++++++
# Helper function
#++++++++++++++++++++++
# Download and parse webpage
html_to_text <- function(url) {
  library(RCurl)
  library(XML)
  # download html
  html.doc <- getURL(url)
  # convert to plain text
  doc <- htmlParse(html.doc, asText = TRUE)
  # "//text()" returns all text outside of HTML tags.
  # We also don’t want text such as style and script codes
  text <- xpathSApply(doc, "//text()[not(ancestor::script)][not(ancestor::style)][not(ancestor::noscript)][not(ancestor::form)]", xmlValue)
  # Format text vector into one character string
  return(paste(text, collapse = " "))
}
# source('http://www.sthda.com/upload/rquery_wordcloud.r')
filePath <- "https://raw.githubusercontent.com/vladwelt/buscador/master/documentos/Romancero%20gitano%20-%20Federico%20Garcia%20Lorca.txt"
res <- rquery.wordcloud(filePath, type = "file", lang = "spanish", max.words = 100)

png("22.png", width = 6, height = 6, units = "in", res = 900)
res <- rquery.wordcloud(filePath, type = "file", lang = "spanish", max.words = 100)
dev.off()
png 
  2

Para un circulo con radio ‘r radio’, el area es, ‘r pi * radio^2’.

Embed a web page

If you have the webshot package [@R-webshot] and PhantomJS installed (see Section @ref(html-widgets)), you can embed any web page in the output document through knitr::include_url(). When you pass a URL of a web page to this function in a code chunk, it will generate an <iframe> (inline frame) if the output format is HTML, and a screenshot of the web page for other output formats. You can view the actual page in the inline frame. For example, Figure @ref(fig:include-url) should show you my homepage if you are reading the online version of this book, otherwise you will see a static screenshot instead.

knitr::include_url('https://yihui.org')

Embed Yihui’s homepage as an iframe or screenshot.

Code chunks and inline R code

An R Markdown document consists of intermingled prose (narratives) and code. There are two types of code in an Rmd document: code chunks and inline R code. Below is a quick example:

radio <- 5.00  # radius of a circle

For a circle with the radius 5, el area es: 49.348022

For a circle with the radius 5, el area es: 49.348022 m2

Formatting

The greatest strength of the Markdown language is that its simplicity makes it very easy to read and write even to newcomers. This is its key design principle, as outlined by the creator of the original Markdown language:

A Markdown-formatted document should be publishable as-is, as plain text, without looking like it’s been marked up with tags or formatting instructions.

— John Gruber

However, this comes at a cost of customization. Many features of typical word processors are not directly available in Markdown, e.g.,

We leave it to you to decide whether such features are worth your effort. To some degree, Markdown reflects the philosophy of Stoicism: the “natural world” consists of plain text, and you should not be controlled by the desire for (visual) pleasure. Anyway, this chapter offers some tips on how you can customize the appearance and styling of elements in an R Markdown document.

If you need a reminder in the basics of the Markdown language, the R Markdown cheatsheet at https://www.rstudio.com/resources/cheatsheets/ provides a good overview of the basic syntax.

colorize = function(x, color){
  if (knitr::is_latex_output()) {
    sprintf("\\textcolor{%s}{%s}", color, x)
  } else if (knitr::is_html_output()) {
    sprintf("<span style='color: %s;'>%s</span>", color, x)
  } else x
}

We can then use the code in an inline R expression `r colorize("some words in red", "red")`, which will create some words in red (you will not see the red color if you are reading this book printed in black and white).

As mentioned in Section @ref(linebreaks), whitespaces are often meaningless in Markdown. Markdown will also ignore spaces used for indentation by default. However, we may want to keep the indentation in certain cases, e.g., in verses and addresses. In these situations, we can use line blocks by starting the line with a vertical bar (|). The line breaks and any leading spaces will be preserved in the output. For example:1



| When dollars appear it's a sign
|   that your code does not quite align  
| Ensure that your math  
|   in xaringan hath  
|   been placed on a single long line



Verbatim code chunks

Unfortunately, we cannot wrap the code chunk in another layer of backticks, but instead we must make the code chunk invalid within the source code by inserting `r ''` in the chunk header. This will be evaluated as an inline expression to an empty string by knitr. For this example, the following “code chunk” in the source document:

https://goo.gl/maps/JN5XDTrcFmNFJYVL7

La suma es: 4

La suma es: 12784.25

https://r-charts.com/es/colores/

El Valor de Mercado por Comparación asciende a: 12784.25 euros

El Valor de Mercado por Comparación asciende a: 12784.25 euros

El Valor de Mercado por Comparación asciende a: 1.255555^{6} euros

El Valor de Mercado por Comparación asciende a: 1.255555^{6} euros 

# install.packages("tinytex")
require(tinytex)

\[\begin{align}f(x_1,\dots,x_n;\mu,\sigma) &= \prod_{i=1}^n \varphi_{\mu,\sigma^2}(x_i)\\ &=\frac1{(\sigma\sqrt{2\pi})^n}\prod_{i=1}^n \exp\biggl(-{1 \over 2} \Bigl({x_i-\mu \over \sigma}\Bigr)^2\biggr), \quad(x_1,\ldots,x_n)\in\mathbb{R}^n. \end{align}\]

library(flextable)

https://ardata-fr.github.io/flextable-book/

set_flextable_defaults(big.mark = " ",
  font.size = 10, theme_fun = theme_vanilla,
  padding.bottom = 6,
  padding.top = 6,
  padding.left = 6,
  padding.right = 6,
background.color = "#EFEFEF")

ft <- flextable(airquality[ sample.int(10),])
ft <- add_header_row(ft,
  colwidths = c(4, 2),
  values = c("Air quality", "Time")
)
ft <- theme_vanilla(ft)
ft <- add_footer_lines(ft, "Daily air quality measurements in New York, May to September 1973.")
ft <- color(ft, part = "footer", color = "#666666")
ft <- set_caption(ft, caption = "New York Air Quality Measurements")
ft
CODIGO INSERTAR VIDEO, haz click para expandir 

<object width="640" height="390">
  <param name="movie"
         value="https://www.youtube.com/embed/OFTPDWBSbQM"></param>
  <param name="allowScriptAccess" value="always"></param>
  <embed src="https://www.youtube.com/embed/OFTPDWBSbQM"
         type="application/x-shockwave-flash"
         allowscriptaccess="always"
         width="640" height="390"></embed>
</object>



Redes geoespaciales ordenadas en R

https://sfnetworks.github.io/useR2021/slides#8

# install.packages("sfnetworks")
library(tidyverse)
library(sfnetworks)
roxel %>%
 as_tibble() %>%
 select(name, type)
## # A tibble: 851 x 2
##    name                  type       
##    <chr>                 <fct>      
##  1 Havixbecker Strasse   residential
##  2 Pienersallee          secondary  
##  3 Schulte-Bernd-Strasse residential
##  4 <NA>                  path       
##  5 Welsingheide          residential
##  6 <NA>                  footway    
##  7 <NA>                  footway    
##  8 <NA>                  path       
##  9 <NA>                  track      
## 10 <NA>                  track      
## # ... with 841 more rows
library(sf)
roxel
Simple feature collection with 851 features and 2 fields
Geometry type: LINESTRING
Dimension:     XY
Bounding box:  xmin: 7.522594 ymin: 51.94151 xmax: 7.546705 ymax: 51.9612
Geodetic CRS:  WGS 84
# A tibble: 851 x 3
   name                  type                                           geometry
 * <chr>                 <fct>                                  <LINESTRING [°]>
 1 Havixbecker Strasse   residential      (7.533722 51.95556, 7.533461 51.95576)
 2 Pienersallee          secondary   (7.532442 51.95422, 7.53236 51.95377, 7.53~
 3 Schulte-Bernd-Strasse residential (7.532709 51.95209, 7.532823 51.95239, 7.5~
 4 <NA>                  path        (7.540063 51.94468, 7.539696 51.94479, 7.5~
 5 Welsingheide          residential       (7.537673 51.9475, 7.537614 51.94562)
 6 <NA>                  footway     (7.543791 51.94733, 7.54369 51.94686, 7.54~
 7 <NA>                  footway           (7.54012 51.94478, 7.539931 51.94514)
 8 <NA>                  path        (7.53822 51.94546, 7.538131 51.94549, 7.53~
 9 <NA>                  track       (7.540063 51.94468, 7.540338 51.94468, 7.5~
10 <NA>                  track       (7.5424 51.94599, 7.54205 51.94629, 7.5419~
# ... with 841 more rows
lines = roxel %>%
  st_transform(3035)
plot(st_geometry(lines))

Cómo Automatizar Scripts de R en Windows

https://anderfernandez.com/blog/automatizar-scripts-de-r-en-windows-y-mac/

La automatización de tareas la realizaremos con el paquete taskscheduleR. Por lo tanto, lo primero será descargar e instalar este paquete.

# install.packages("taskscheduleR")
library(taskscheduleR)
# taskscheduleR:::taskschedulerAddin()
library(lubridate)
ahora <- Sys.time()
ahora <- gsub(" |:","_",ahora)
write.csv(mtcars,paste0(ahora,"mtcars.csv"))
# ABRIR ADDINS

# library(taskscheduleR)
# fichero <- "E:\\Automatizar_mtcars.R"
# 
# taskscheduler_create(taskname = "mtcars", 
#                      rscript = fichero,
#                      schedule = "MINUTE", 
#                      starttime = format(Sys.time(), "%H:%M:%S"), 
#                      startdate = format(Sys.time(), "%d/%m/%Y"))

INTRODUCIR LAS LÍNEAS DE UN SCRIPT EN R MARKDOWN SIN EJECUTARLO # Función de leer e imprimir un fichero .R

script <- readLines("RstudioBAS.R") 
cat(script, sep = '\n') 
  # **************************************
  # Manejo de archivos y Workspace en R
  # **************************************
  
  # Definir el directorio de trabajo (working directory)
  setwd("c:/")
  
  # Obtener lista completa de archivos en el Working directory
  dir()
  
  # Obtener lista completa de data frames R en el Working directory
  dir(pattern=".Rda")
  
  # Obtener lista completa de scripts R en el Working directory
  # El dolar $ implica que es final exacto de la extensión,
  # no presenta archivos Rda por ejemplo
  dir(pattern=".R$")
  
  # Leer un archivo del working directory: no se pone path, pues
  # es el path del setwd
  
  load("saheartbis.Rda")
  
  # Guardar  un archivo en el working directory: no se pone path, pues
  # es el path del setwd
  
  save(saheartbis,file="saheartbis.Rda")
  
  # Grabar todos los objetos del Global Environment (Workspace)  
  # en un archivo
  
  save.image(file="todosobjetos.RData")
  
  # Leer  objetos del Workspace de un archivo
  
  load("todosobjetos.RData")
  
  # Borrar (remover) objetos que ocupan RAM en el Global Environment
  
  rm(saheartbis)
  
  # Borrar todos los objetos del Global Environment 
  # y borrar basura de la memoria
  
  rm(list=ls()) 
  gc()
  
  # Borrar todos los plots
  graphics.off()
  
  # También vale 
  # dev.off()

Como se puede ver en el código anterior…

Here is the first Div.

# {r fig.height = 3, fig.width = 5}
plot (pressure)

And this block will be put on the right:

# {r fig.height = 3, fig.width = 5, fig.align = "center"}
plot (pressure)


Here is the first Div.

# {r fig.width = 5, fig.asp = 0.62}
plot (pressure)

And this block will be put on the right:

#{r out.width = "40%"}
plot (pressure)


insertar el scrip del INE directamente en la WEB


aside element



PARA VER EL CODIGO, haz click para expandir 



<style>
aside {

width : 30%;
padding-left : 15px;
margin-left : 15px;
float : right ;
font-style : italic;
background-color : lightgray;
}
</style>

<body>

<h5>The aside element - Styled with CSS</h5>
<p>My family and I visited The Epcot center this
summer. The weather was nice, and Epcot was amazing!
I had a great summer together with my family! </p>

<aside>
<p>The Epcot center is a theme park at Walt Disney
World Resort featuring exciting attractions,
international pavilions, award-winning fireworks and
seasonal special events. </p>
</aside>

<p>My otra cantidad de palabras para probar que
se puede hacer con este bloque de codigo</p>
</body>



The aside element - Styled with CSS

My family and I visited The Epcot center this summer. The weather was nice, and Epcot was amazing! I had a great summer together with my family!

My otra cantidad de palabras para probar que se puede hacer con este bloque de codigo








tipos de lineas para insertar rectas en graficos https://psyteachr.github.io/introdataviz/additional-customisation-options.html

linetype = 0 linetype = “blank”
linetype = 1 linetype = “solid”
linetype =2 linetype =“dashed”
linetype = 3 linetype = “dotted”
linetype = 4 linetype = “dotdash”
linetype = 5 linetype = “longdash”
linetype = 6 linetype = “twodash”



# install.packages("parchwork")
library(patchwork)
library(tidyverse)
load("C:/Users/polo/OneDrive/ZZZ/VisualizacionDatosR/dat_long.RData")

Here is the first Div.

ggplot(dat_long, aes(x = acc)) +
  geom_histogram(binwidth = 1, fill = "white", color = "black") +
  scale_x_continuous(name = "Accuracy (0-100)") +
  geom_vline(xintercept = 80, linetype = 2, color = "red", size = 1.5)

And this block will be put on the right:

ggplot(dat_long, aes(x = condition, y = acc)) +
  geom_boxplot() +
  geom_hline(yintercept = 80, linetype = 3, color = "blue", size = 2)



a <- ggplot(dat_long, aes(x = rt)) +
  geom_histogram(binwidth = 10, fill = "white", color = "black") +
  scale_x_continuous(name = "Reaction time (ms)")

b <- ggplot(dat_long, aes(x = rt)) +
  geom_density()+
  scale_x_continuous(name = "Reaction time (ms)")

a + inset_element(b, left = 0.6, bottom = 0.6, right = 1, top = 1)



# install.packages("leafem")
library(leaflet)
library(leaflegend)
library(leafem)

img <- "https://static.wixstatic.com/media/e5cab7_6b500fcdb224460da6678532531725b4~mv2.png/v1/fill/w_176,h_155,al_c,q_85,usm_0.66_1.00_0.01/logo_aeviu.webp"


# leaflet() %>% addTiles() %>% addLogo(img, url = "https://www.r-project.org/logo/")
leaflet() %>% addTiles() %>% addLogo(img, url = "https://www.aeviu.es/")



https://github.com/coolbutuseless/ggpattern

# install.packages("remotes")
# remotes::install_github("coolbutuseless/ggpattern")
library (remotes)
library(ggpattern)

Tome una parcela existente que contenga una geometría con un área rellenable, por ejemplo.geom_col() Utilice la versión de la geometría, por ejemplo. En lugar de {ggpattern}ggpattern::geom_col_pattern()ggplot2::geom_col() Establezca la estética a su elección de patrón, por ejemplo, y establezca otras opciones utilizando la estéticapatternpattern = ’stripe’pattern_*

df <- data.frame(level = c("a", "b", "c", 'd'), outcome = c(2.3, 1.9, 3.2, 1))

ggplot(df) +
  geom_col_pattern(
    aes(level, outcome, pattern_fill = level), 
    pattern = 'stripe',
    fill    = 'white',
    colour  = 'black'
  ) +
  theme_bw(18) +
  theme(legend.position = 'none') + 
  labs(
    title    = "ggpattern::geom_col_pattern()",
    subtitle = "pattern = 'stripe'"
  ) +
  coord_fixed(ratio = 1/2)



Parámetros para patrones basados en geometría https://coolbutuseless.github.io/package/ggpattern/articles/geometry-based-pattern-parameters.html https://coolbutuseless.github.io/package/ggpattern/articles/pattern-image.html https://coolbutuseless.github.io/package/ggpattern/articles/pattern-stripe.html https://coolbutuseless.github.io/package/ggpattern/articles/pattern-magick.html https://coolbutuseless.github.io/package/ggpattern/articles/pattern-placeholder.html https://coolbutuseless.github.io/package/ggpattern/articles/developing-patterns.html https://coolbutuseless.github.io/package/ggpattern/articles/gganimate.html



suppressPackageStartupMessages({
  library(dplyr)
  library(ggplot2)
  library(ggpattern)
})

#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Simple testing data
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
df1 <- data.frame(
  trt     = c("a", "b", "c"), 
  outcome = c(2.3, 1.9, 3.2),
  stringsAsFactors = FALSE
)

ggplot(df1, aes(trt, outcome)) +
  geom_col_pattern(
    aes(fill = trt), 
    pattern      = 'placeholder', 
    pattern_type = 'bear',
    colour       = 'black'
  ) +
  theme_bw(15) +
  labs(
    title    = "ggpattern::geom_col_pattern()",
    subtitle = "pattern='placeholder', pattern_type='bear'"
  ) +
  theme(legend.position = 'none') +
  coord_fixed(ratio = 1/2)

https://rdrr.io/github/rOpenSpain/caRtociudad/man/get_cartociudad_area.html

polo1<-c(
  "40,5530339, -4,016714",
  "40,5421333, -4,0207559",
  "40,5480905, -4,0059137",
  "40,5425003, -4,0211858",
  "40,5500285, -4,0041974",
  "40,5416315, -4,0215055")
library(readxl)
library(data.table)         # Manipulacion
library(knitr)              # RMardown
library(rJava)
library(xlsx)
library(xlsxjars)

tablas<-read_excel("C:/Users/polo/Desktop/AEVIU/0002/0002piso.xlsm", sheet="elemento", range=("A63:B69"))

tablas<-as.data.frame(tablas)
# tablas<-as.vector(tablas)
tablas
       lat       lng
1 40.55303 -4.016714
2 40.54920 -4.008311
3 40.54809 -4.005914
4 40.54250 -4.021186
5 40.55003 -4.004197
6 40.54163 -4.021505
# coord <- data.frame(lat = c(40.57455), lng = c(-4.00189))
# coord[2,] <- c(40.60109, -3.98937)
# coord[3,] <- c(40.60890, -3.94460)
# coord[4,] <- c(40.58923, -3.92543)
# coord[5,] <- c(40.58063, -3.95083)
# coord[6,] <- c(40.58413, -3.99883)
# coord[7,] <- c(40.61589, -4.00917)
# coord
#pegar como vector vertical

# direcciones<-c("calle Prado Ibarra 23, 28270",
#    "calle Sacramientos  9 C, 28270",
#    "Travesia Viñas Viejas 7, 28270",
#    "calle Calderon de la Barca 29, 28270",
#    "calle Viñas Viejas 45, 28270",
#    "calle Calderon de la Barca 47, 28270")
#  
# direcciones[1]
# direcciones[2]
# direcciones[3]
# direcciones[4]
# direcciones[5]
# direcciones[6]

standard usage res <- cartociudad_geocode(full_address = “plaza de cascorro 11, 28005 madrid”)

km 41 of A-23 motorway res <- cartociudad_geocode(“A-23 41”)

specific usage (see References for details) res <- cartociudad_geocode(“A-23 41”, type = “portal”, id = “600000000045”, portal = 41)

vectorized call addresses <- paste(“A-23”, 1:10) res <- lapply(addresses, cartociudad_geocode, on.error = “warn”)

# copiando de la hoja de calculo #########################################################################

direcciones<-read_excel("C:/Users/polo/Desktop/AEVIU/0002/0002piso.xlsm", sheet="elemento", range=("c63:d69"))
direcciones
# A tibble: 6 x 2
  tes                        cp
  <chr>                   <dbl>
1 Prado Ibarra 23         28270
2 Gregorio Panadero 28    28270
3 Viñas Viejas 7          28270
4 Calderon de la Barca 29 28270
5 Viñas Viejas 45         28270
6 Calderon de la Barca 47 28270
dir1<-paste(direcciones$tes[1], direcciones$cp[1], sep=", ")
dir2<-paste(direcciones$tes[2], direcciones$cp[2], sep=", ")
dir3<-paste(direcciones$tes[3], direcciones$cp[3], sep=", ")
dir4<-paste(direcciones$tes[4], direcciones$cp[4], sep=", ")
dir5<-paste(direcciones$tes[5], direcciones$cp[5], sep=", ")
dir6<-paste(direcciones$tes[6], direcciones$cp[6], sep=", ")

dir1
[1] "Prado Ibarra 23, 28270"
dir2
[1] "Gregorio Panadero 28, 28270"
dir3
[1] "Viñas Viejas 7, 28270"
dir4
[1] "Calderon de la Barca 29, 28270"
dir5
[1] "Viñas Viejas 45, 28270"
dir6
[1] "Calderon de la Barca 47, 28270"
library(caRtociudad)
library(ggmap)
library(doParallel)
library(leaflet)
library(remotes)
# library(rOpenSpain)
# install.packages("remotes")
# remotes::install_github("rOpenSpain/caRtociudad")

my.address1 <- cartociudad_geocode(dir1)
my.address2 <- cartociudad_geocode(dir2)
my.address3 <- cartociudad_geocode(dir3)
my.address4 <- cartociudad_geocode(dir4)
my.address5 <- cartociudad_geocode(dir5)
my.address6 <- cartociudad_geocode(dir6)

cuadro<-rbind(my.address1, my.address2, my.address3, my.address4, my.address5, my.address6)

cuadro<-as.matrix(cuadro)
#cuadro<-t(cuadro)
cuadro
     id              province comunidadAutonoma     muni          type    
[1,] "280440000028"  "Madrid" "Comunidad de Madrid" "Colmenarejo" "portal"
[2,] "280440000103"  "Madrid" "Comunidad de Madrid" "Colmenarejo" "portal"
[3,] "280440000111"  "Madrid" "Comunidad de Madrid" "Colmenarejo" "portal"
[4,] "2280440670409" "Madrid" "Comunidad de Madrid" "Colmenarejo" "portal"
[5,] "280440000111"  "Madrid" "Comunidad de Madrid" "Colmenarejo" "portal"
[6,] "2280440681488" "Madrid" "Comunidad de Madrid" "Colmenarejo" "portal"
     address              postalCode poblacion    
[1,] "PRADO IBARRA (DEL)" "28270"    "Colmenarejo"
[2,] "GREGORIO PANADERO"  "28270"    "Colmenarejo"
[3,] "VIÑAS VIEJAS"       "28270"    "Colmenarejo"
[4,] "CALDERON BARCA"     "28270"    "Colmenarejo"
[5,] "VIÑAS VIEJAS"       "28270"    "Colmenarejo"
[6,] "CALDERON BARCA"     "28270"    "Colmenarejo"
     geom                                           tip_via lat       
[1,] "POINT (-4.01671406657353 40.55303398891171)"  "CALLE" "40.55303"
[2,] "POINT (-4.008147011475415 40.54918279414671)" "CALLE" "40.54918"
[3,] "POINT (-4.017183557675118 40.54034451428396)" "CALLE" "40.54034"
[4,] "POINT (-4.005854696900138 40.54806652816524)" "Calle" "40.54807"
[5,] "POINT (-4.020174038201369 40.5411915947124)"  "CALLE" "40.54119"
[6,] "POINT (-4.007418841794464 40.54759463145535)" "Calle" "40.54759"
     lng         portalNumber stateMsg                                 state
[1,] "-4.016714" "23"         "Resultado exacto de la búsqueda"        "1"  
[2,] "-4.008147" "28"         "Resultado exacto de la búsqueda"        "1"  
[3,] "-4.017184" "7"          "Resultado exacto de la búsqueda"        "1"  
[4,] "-4.005855" "29"         "Resultado exacto de la búsqueda"        "1"  
[5,] "-4.020174" "43"         "Portal no encontrado. Par más cercano." "2"  
[6,] "-4.007419" "47"         "Resultado exacto de la búsqueda"        "1"  
     countryCode
[1,] "011"      
[2,] "011"      
[3,] "011"      
[4,] "011"      
[5,] "011"      
[6,] "011"      
my.address1$lng
[1] -4.016714
# DE CADA LONGITUD-LATITUD SE OBTINEN TOTODOS LOS DATOS INCLUIDOLA REFERECIA CATASTRAL Y WEB

tesp10<-cartociudad_get_location_info(my.address1$lat, my.address1$lng, 2020, info.source = "cadastre")
tesp11<-cartociudad_get_location_info(my.address1$lat, my.address1$lng, 2020, info.source = "reverse")
tesp10<-t(tesp10)
tesp11<-t(tesp11)
# tesp10
# tesp11
tesp1<-cbind(tesp10, tesp11)
tesp1
     ref.catastral   
[1,] "4000538VK1940S"
     url.ref.catastral                                                                                                
[1,] "https://www1.sedecatastro.gob.es/CYCBienInmueble/OVCListaBienes.aspx?del=28&muni=44&rc1=4000538&rc2" [truncated]
     tipo tipo.via nombre.via           num.via num.via.id     municipio    
[1,] NULL "CALLE"  "PRADO IBARRA (DEL)" 23      "280440127635" "Colmenarejo"
     provincia cod.postal
[1,] "Madrid"  "28270"   
tesp20<-cartociudad_get_location_info(my.address2$lat, my.address2$lng, 2020, info.source = "cadastre")
tesp21<-cartociudad_get_location_info(my.address2$lat, my.address2$lng, 2020, info.source = "reverse")
tesp20<-t(tesp20)
tesp21<-t(tesp21)
#tesp20
#tesp21
tesp2<-cbind(tesp20, tesp21)
tesp2
     tipo tipo.via nombre.via          num.via num.via.id     municipio    
[1,] NULL "CALLE"  "GREGORIO PANADERO" 28      "280440674508" "Colmenarejo"
     provincia cod.postal
[1,] "Madrid"  "28270"   
tesp30<-cartociudad_get_location_info(my.address3$lat, my.address3$lng, 2020, info.source = "cadastre")
tesp31<-cartociudad_get_location_info(my.address3$lat, my.address3$lng, 2020, info.source = "reverse")
tesp30<-t(tesp30)
tesp31<-t(tesp31)
# tesp30
# tesp31
tesp3<-cbind(tesp30, tesp31)
tesp3
     tipo tipo.via nombre.via     num.via num.via.id     municipio    
[1,] NULL "CALLE"  "VIÑAS VIEJAS" 7       "280440674427" "Colmenarejo"
     provincia cod.postal
[1,] "Madrid"  "28270"   
tesp40<-cartociudad_get_location_info(my.address4$lat, my.address4$lng, 2020, info.source = "cadastre")
tesp41<-cartociudad_get_location_info(my.address4$lat, my.address4$lng, 2020, info.source = "reverse")
tesp40<-t(tesp40)
tesp41<-t(tesp41)
# tesp40
# tesp41
tesp4<-cbind(tesp40, tesp41)
tesp4
     tipo     tipo.via nombre.via       num.via num.via.id      municipio    
[1,] "portal" "Calle"  "CALDERON BARCA" 29      "2280440670409" "Colmenarejo"
     provincia cod.postal
[1,] "Madrid"  "28270"   
tesp50<-cartociudad_get_location_info(my.address4$lat, my.address4$lng, 2020, info.source = "cadastre")
tesp51<-cartociudad_get_location_info(my.address5$lat, my.address5$lng, 2020, info.source = "reverse")
tesp50<-t(tesp50)
tesp51<-t(tesp51)
# tesp50
# tesp51
tesp5<-cbind(tesp50, tesp51)
tesp5
     tipo tipo.via nombre.via     num.via num.via.id     municipio    
[1,] NULL "CALLE"  "VIÑAS VIEJAS" 43      "280440674688" "Colmenarejo"
     provincia cod.postal
[1,] "Madrid"  "28270"   
tesp60<-cartociudad_get_location_info(my.address6$lat, my.address6$lng, 2020, info.source = "cadastre")
tesp61<-cartociudad_get_location_info(my.address6$lat, my.address6$lng, 2020, info.source = "reverse")
tesp60<-t(tesp60)
tesp61<-t(tesp61)
# tesp60
# tesp61
tesp6<-cbind(tesp60, tesp61)
tesp6
     ref.catastral   
[1,] "4993698VK1859S"
     url.ref.catastral                                                                                                
[1,] "https://www1.sedecatastro.gob.es/CYCBienInmueble/OVCListaBienes.aspx?del=28&muni=44&rc1=4993698&rc2" [truncated]
     tipo     tipo.via nombre.via       num.via num.via.id      municipio    
[1,] "portal" "Calle"  "CALDERON BARCA" 47      "2280440681488" "Colmenarejo"
     provincia cod.postal
[1,] "Madrid"  "28270"   
testigosparcelas<-rbind(tesp1, tesp6)
testigosparcelas
     ref.catastral   
[1,] "4000538VK1940S"
[2,] "4993698VK1859S"
     url.ref.catastral                                                                                                
[1,] "https://www1.sedecatastro.gob.es/CYCBienInmueble/OVCListaBienes.aspx?del=28&muni=44&rc1=4000538&rc2" [truncated]
[2,] "https://www1.sedecatastro.gob.es/CYCBienInmueble/OVCListaBienes.aspx?del=28&muni=44&rc1=4993698&rc2" [truncated]
     tipo     tipo.via nombre.via           num.via num.via.id     
[1,] NULL     "CALLE"  "PRADO IBARRA (DEL)" 23      "280440127635" 
[2,] "portal" "Calle"  "CALDERON BARCA"     47      "2280440681488"
     municipio     provincia cod.postal
[1,] "Colmenarejo" "Madrid"  "28270"   
[2,] "Colmenarejo" "Madrid"  "28270"   
# area con un radio en metros

get_cartociudad_area(40.4873817, -3.3826135, 500)
   longitude latitude
1    -3.3860  40.4849
2    -3.3868  40.4850
3    -3.3877  40.4856
4    -3.3878  40.4866
5    -3.3864  40.4907
6    -3.3863  40.4908
7    -3.3811  40.4916
8    -3.3810  40.4915
9    -3.3789  40.4894
10   -3.3784  40.4855
11   -3.3786  40.4854
12   -3.3860  40.4849
res <- cartociudad_get_route(c(39.48,-0.37),
   c(39.484336,-0.358171),
   vehicle = "car")
res
$bbox
[1] -0.3711818 39.4799157 -0.3582531 39.4853696

$distance
[1] 1764.33

$found
[1] TRUE

$from
[1] 39.48 -0.37

$geom
[1] "e}}oFnggA?NVzAQEaAk@E@a@rAUOa@lB[SeFmB_EcByD}AkAk@uB_A@eAJsDUY@U@UBmALgE@S?YPaFEUMIMGWQKGIEWE^eBLe@WKAWbA{EaBo@nAqGDWbAmFDWx@eEMC"

$info
$info$routeFound
[1] TRUE

$info$took
[1] 24

$info$tookGeocoding
[1] 0


$instructionsData
        bbox1    bbox2      bbox3    bbox4
1  -0.3700731 39.48003 -0.3699971 39.48003
2  -0.3707131 39.47992 -0.3702801 39.48065
3  -0.3711818 39.48082 -0.3693424 39.48438
4  -0.3690265 39.48480 -0.3645113 39.48500
5  -0.3644751 39.48496 -0.3644751 39.48496
6  -0.3643891 39.48508 -0.3634051 39.48531
7  -0.3623087 39.48488 -0.3623087 39.48488
8  -0.3620671 39.48456 -0.3592664 39.48537
9  -0.3582707 39.48427 -0.3582707 39.48427
10 -0.3582531 39.48435 -0.3582531 39.48435
                                               description      dest1    dest2
1                   Continúe por CALLE VUELTA DEL RUISEÑOR -0.3700731 39.48003
2                  Gire justo a la derecha por CALLE FLORA -0.3706364 39.48065
3               Gire justo a la derecha por CALLE ALBORAYA -0.3693424 39.48438
4                     Gire a la derecha por CALLE MOLINELL -0.3645113 39.48489
5                                Continúe por PRIMADO REIG -0.3644751 39.48496
6               Continúe por CALLE DOCTOR VICENTE ZARAGOZA -0.3634051 39.48521
7                Gire a la izquierda por CALLE EMILIO BARO -0.3623087 39.48488
8  Gire a la derecha por CALLE REVERENDO RAFAEL TRAMOYERES -0.3592664 39.48456
9                    Gire a la izquierda por CALLE MISTRAL -0.3582707 39.48427
10                                        Objetivo logrado -0.3582531 39.48435
   distance indication      orig1    orig2
1        48          0 -0.3699971 39.48003
2       160          3 -0.3700731 39.48003
3       496          3 -0.3706364 39.48065
4       408          2 -0.3693424 39.48438
5        15          0 -0.3645113 39.48489
6       219          0 -0.3644751 39.48496
7        58         -2 -0.3634051 39.48521
8       348          2 -0.3623087 39.48488
9         8         -2 -0.3592664 39.48456
10        0          4 -0.3582707 39.48427

$time
[1] 211706

$to
[1] 39.484336 -0.358171
# install.packages("ggmap")
library(ggmap)

soria <- cartociudad_geocode("plaza de san esteban, soria")
soria_map <- cartociudad_get_map(c(soria$lat, soria$lng), 1)
ggmap::ggmap(soria_map)

info<-cartociudad_get_location_info(40.4873817, -3.3826135, 2021, info.source = "cadastre")
info
$ref.catastral
[1] "7619501VK6872S"

$url.ref.catastral
[1] "https://www1.sedecatastro.gob.es/CYCBienInmueble/OVCListaBienes.aspx?del=28&muni=5&rc1=7619501&rc2=VK6872S"



library(leaflet)
library(mapSpain)

catastromap <-
  leaflet(height=450, width=800) %>%
  setView(
    lat = 40.5594660,
    lng = -3.9766931,
    zoom = 15
  ) %>%
  
  addProviderEspTiles(provider = "IGNBase.TodoNoFondo")%>%
  addProviderEspTiles(provider = "Catastro.CadastralParcel")%>%
  # addProviderEspTiles(provider = "Catastro.CadastralZoning")%>%
  # addProviderEspTiles(provider = "Catastro.Address")
  # addProviderEspTiles(provider = "Catastro.Parcela")
  addProviderEspTiles(provider = "Catastro.Building")
  # addProviderEspTiles(provider = "Cartociudad.CodigosPostales")

catastromap
library(sf)
library(leaflet)
library(dplyr)
library(leafem)
library(htmltools)
library(mapSpain)
library(data.table)  
# install.packages("xlsx")
library(xlsx)


mapatabla<-read.xlsx("C:/Users/polo/OneDrive/ZZZ/vivienda/tablasNAtotal.xlsx", sheetName="Sheet1", rowIndex=c(1, 2:610))

mapatabla<-as.data.frame(mapatabla)
# View(mapatabla)

leaflet(mapatabla, height=450, width=800, options = leafletOptions(minZoom = 5, maxZoom = 20)) %>% addTiles(group = "OSM") %>% addTiles() %>%
addMarkers(~long, ~lat, popup = ~htmlEscape(direccion),
           clusterOptions = markerClusterOptions(maxClusterRadius = 5))%>%
addProviderEspTiles(provider = "IGNBase.TodoNoFondo")%>%
addProviderEspTiles(provider = "Catastro.CadastralParcel")%>%
  # addProviderEspTiles(provider = "Catastro.CadastralZoning")%>%
  # addProviderEspTiles(provider = "Catastro.Address")
  # addProviderEspTiles(provider = "Catastro.Parcela")
addProviderEspTiles(provider = "Catastro.Building")%>%
  # addProviderEspTiles(provider = "Cartociudad.CodigosPostales")
addMeasure()
#   https://github.com/Leaflet/Leaflet.markercluster#customising-the-clustered-markers

#   m %>% addMarkers(popup= ~paste("Hola mundo", size) )



# install.packages("E:/curl-7.81.0.tar.gz", repos = NULL, type = "source")
library(curl)
## Using libcurl 7.64.1 with Schannel
## 
## Attaching package: 'curl'
## The following object is masked from 'package:readr':
## 
##     parse_date
req <- curl_fetch_memory("https://eu.httpbin.org/get?foo=123")
str(req)
## List of 7
##  $ url        : chr "https://eu.httpbin.org/get?foo=123"
##  $ status_code: int 200
##  $ type       : chr "application/json"
##  $ headers    : raw [1:230] 48 54 54 50 ...
##  $ modified   : POSIXct[1:1], format: NA
##  $ times      : Named num [1:6] 0 0.0283 0.1286 0.3906 0.5292 ...
##   ..- attr(*, "names")= chr [1:6] "redirect" "namelookup" "connect" "pretransfer" ...
##  $ content    : raw [1:361] 7b 0a 20 20 ...

http://madrid.r-es.org/wp-content/uploads/2016/05/caRtociudad_20160512.pdf

https://github.com/thinkr-open/remedy

# remotes::install_github("ThinkR-open/remedy")
library(remedy)

https://www.vishalkatti.com/posts/2021-07-17-programmingwithdplyr/

Utilizaremos los datos de Ventas de Vivienda de Texas, disponibles como un tibble en el paquete popular como datos de referencia. Contiene información mensual sobre el mercado de la vivienda en Texas proporcionada por el centro de bienes raíces TAMU, https://www.recenter.tamu.edu/. Tiene 8602 observaciones y 9 variables.ggplot2

txhousing <- ggplot2::txhousing
dplyr::glimpse(txhousing)
## Rows: 8,602
## Columns: 9
## $ city      <chr> "Abilene", "Abilene", "Abilene", "Abilene", "Abilene", "Abil~
## $ year      <int> 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, ~
## $ month     <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 1, 2, 3, 4, 5, 6, 7, ~
## $ sales     <dbl> 72, 98, 130, 98, 141, 156, 152, 131, 104, 101, 100, 92, 75, ~
## $ volume    <dbl> 5380000, 6505000, 9285000, 9730000, 10590000, 13910000, 1263~
## $ median    <dbl> 71400, 58700, 58100, 68600, 67300, 66900, 73500, 75000, 6450~
## $ listings  <dbl> 701, 746, 784, 785, 794, 780, 742, 765, 771, 764, 721, 658, ~
## $ inventory <dbl> 6.3, 6.6, 6.8, 6.9, 6.8, 6.6, 6.2, 6.4, 6.5, 6.6, 6.2, 5.7, ~
## $ date      <dbl> 2000.000, 2000.083, 2000.167, 2000.250, 2000.333, 2000.417, ~
select_raw <- function(df, var) {
  dplyr::select(.data = df, {{var}}) %>%     # embrace of curly-curly {{}} brackets
    head()                                   # para limitar el número de filas de salida en este ejemplo
}
select_raw(txhousing, sales)                 # pasar un solo nombre sin procesar
## # A tibble: 6 x 1
##   sales
##   <dbl>
## 1    72
## 2    98
## 3   130
## 4    98
## 5   141
## 6   156
select_raw(txhousing, c(sales, volume))      # pasar un vector de nombres sin procesar para múltiples columnas
## # A tibble: 6 x 2
##   sales   volume
##   <dbl>    <dbl>
## 1    72  5380000
## 2    98  6505000
## 3   130  9285000
## 4    98  9730000
## 5   141 10590000
## 6   156 13910000

En este método, pasamos la condición como una expresión cruda / desnuda.sales > 8000

filter_raw <- function(df, cond) {
  dplyr::filter(.data = df, {{cond}})        # embrace of curly-curly {{}} brackets
}

filter_raw(txhousing, sales > 8000)          # Pass a single raw criterion
## # A tibble: 10 x 9
##    city     year month sales     volume median listings inventory  date
##    <chr>   <int> <int> <dbl>      <dbl>  <dbl>    <dbl>     <dbl> <dbl>
##  1 Houston  2006     5  8040 1602621368 151200    35398       5.5 2006.
##  2 Houston  2006     6  8628 1795898108 155200    36281       5.6 2006.
##  3 Houston  2013     5  8439 2121508529 186100    20526       3.3 2013.
##  4 Houston  2013     7  8468 2168720825 187800    21497       3.3 2014.
##  5 Houston  2013     8  8155 2083377894 186700    21366       3.3 2014.
##  6 Houston  2014     6  8391 2342443127 211200    19725       2.9 2014.
##  7 Houston  2014     7  8391 2278932511 199700    20214       3   2014.
##  8 Houston  2014     8  8167 2195184825 202400    20007       2.9 2015.
##  9 Houston  2015     6  8449 2490238594 222400    22311       3.2 2015.
## 10 Houston  2015     7  8945 2568156780 217600    23875       3.4 2016.

Pasar múltiples criterios sin procesar usando … argumento Para pasar múltiples criterios sin procesar, podemos usar el argumento….

my_filter <- function(df, ...) { 
  dplyr::filter(.data = df, ...)                # pass the dots argument
  }

my_filter(txhousing, sales > 8000, year > 2010) # pasar múltiples criterios sin procesar
## # A tibble: 8 x 9
##   city     year month sales     volume median listings inventory  date
##   <chr>   <int> <int> <dbl>      <dbl>  <dbl>    <dbl>     <dbl> <dbl>
## 1 Houston  2013     5  8439 2121508529 186100    20526       3.3 2013.
## 2 Houston  2013     7  8468 2168720825 187800    21497       3.3 2014.
## 3 Houston  2013     8  8155 2083377894 186700    21366       3.3 2014.
## 4 Houston  2014     6  8391 2342443127 211200    19725       2.9 2014.
## 5 Houston  2014     7  8391 2278932511 199700    20214       3   2014.
## 6 Houston  2014     8  8167 2195184825 202400    20007       2.9 2015.
## 7 Houston  2015     6  8449 2490238594 222400    22311       3.2 2015.
## 8 Houston  2015     7  8945 2568156780 217600    23875       3.4 2016.

Pasar criterios individuales como una cadena de caracteres De forma predeterminada, no acepta condiciones como cadenas de caracteres. A continuación se muestra un ejemplo que resulta en un errordplyr::filter()

my_filter_string <- function(df, cond) {
  dplyr::filter(.data = df, eval(parse(text = cond)))   # convert text to raw criterion
}

my_filter_string(txhousing, "sales > 8000")             # pass single text string as criteria
## # A tibble: 10 x 9
##    city     year month sales     volume median listings inventory  date
##    <chr>   <int> <int> <dbl>      <dbl>  <dbl>    <dbl>     <dbl> <dbl>
##  1 Houston  2006     5  8040 1602621368 151200    35398       5.5 2006.
##  2 Houston  2006     6  8628 1795898108 155200    36281       5.6 2006.
##  3 Houston  2013     5  8439 2121508529 186100    20526       3.3 2013.
##  4 Houston  2013     7  8468 2168720825 187800    21497       3.3 2014.
##  5 Houston  2013     8  8155 2083377894 186700    21366       3.3 2014.
##  6 Houston  2014     6  8391 2342443127 211200    19725       2.9 2014.
##  7 Houston  2014     7  8391 2278932511 199700    20214       3   2014.
##  8 Houston  2014     8  8167 2195184825 202400    20007       2.9 2015.
##  9 Houston  2015     6  8449 2490238594 222400    22311       3.2 2015.
## 10 Houston  2015     7  8945 2568156780 217600    23875       3.4 2016.

Pasar varios criterios como vector de caracteres ¿Qué pasa si desea pasar varios criterios como un vector de cadena? En tal situación, debemos combinar todas las condiciones de cadena en una sola condición de cadena larga usando . Combina todos los criterios en un solo criterio largo, pero sigue siendo una cadena de texto.paste0(…, collapse = " & “)paste0(”(“, cond,”)“, collapse =” & ")

my_filter_strings <- function(df, cond) { 
  filter_text <- paste0("(", cond, ")", collapse = " & ")   # combine all criteria
  message("Filter Condition: ", filter_text)                # (OPTIONAL) show the combined filter string
  dplyr::filter(.data = df, eval(parse(text = filter_text)))# convert text to raw criterion
  }

my_filter_criteria <- c("sales > 8000", "year > 2010")
my_filter_strings(txhousing, my_filter_criteria)
## Filter Condition: (sales > 8000) & (year > 2010)
## # A tibble: 8 x 9
##   city     year month sales     volume median listings inventory  date
##   <chr>   <int> <int> <dbl>      <dbl>  <dbl>    <dbl>     <dbl> <dbl>
## 1 Houston  2013     5  8439 2121508529 186100    20526       3.3 2013.
## 2 Houston  2013     7  8468 2168720825 187800    21497       3.3 2014.
## 3 Houston  2013     8  8155 2083377894 186700    21366       3.3 2014.
## 4 Houston  2014     6  8391 2342443127 211200    19725       2.9 2014.
## 5 Houston  2014     7  8391 2278932511 199700    20214       3   2014.
## 6 Houston  2014     8  8167 2195184825 202400    20007       2.9 2015.
## 7 Houston  2015     6  8449 2490238594 222400    22311       3.2 2015.
## 8 Houston  2015     7  8945 2568156780 217600    23875       3.4 2016.
my_filter_criteria_with_OR <- c("sales > 8000 | sales < 50", "year > 2010")

# NOTE: OR criteria must be a single string separated by pipe '|' as in example below.
my_filter_strings(txhousing, my_filter_criteria_with_OR)
## Filter Condition: (sales > 8000 | sales < 50) & (year > 2010)
## # A tibble: 315 x 9
##    city         year month sales  volume median listings inventory  date
##    <chr>       <int> <int> <dbl>   <dbl>  <dbl>    <dbl>     <dbl> <dbl>
##  1 Brownsville  2011     1    48 4974408  83300      784      12.6 2011 
##  2 Brownsville  2011     2    47 5558575 101400      776      12.7 2011.
##  3 Brownsville  2011     7    47 4807019  91200      749      13.1 2012.
##  4 Brownsville  2011    12    39 4203440  86800      726      12.4 2012.
##  5 Brownsville  2012     1    43 3892348  85000      791      13.6 2012 
##  6 Brownsville  2012     3    27 2976148  93800      734      13.3 2012.
##  7 Brownsville  2012    11    41 5115393  99000      807      14   2013.
##  8 Brownsville  2013    11    38 4824930 108000      859      13.4 2014.
##  9 Brownsville  2015     1    41 5400796  97000      733      10.7 2015 
## 10 Galveston    2011     1    43 8882961 170000     1015      13.7 2011 
## # ... with 305 more rows

mutate() le permite agregar nuevas columnas o modificar columnas existentes. En el siguiente ejemplo, crearemos una nueva columna a partir de la columna existente. Los nombres de ambas columnas se pueden pasar a la función como nombres sin formato o cadenas de caracteres.volume_in_millionsvolume

mutate_raw <- function(df, new_col_raw, old_col_raw, num = 1) { 
  dplyr::mutate(.data = df, {{new_col_raw}} := {{old_col_raw}}/num) %>% 
    head()
}

txhousing %>% 
  select(city, year, month, volume) %>% 
  mutate_raw(vol_in_millions, volume, 1E6) # pass raw column names w/o quotes
## # A tibble: 6 x 5
##   city     year month   volume vol_in_millions
##   <chr>   <int> <int>    <dbl>           <dbl>
## 1 Abilene  2000     1  5380000            5.38
## 2 Abilene  2000     2  6505000            6.50
## 3 Abilene  2000     3  9285000            9.28
## 4 Abilene  2000     4  9730000            9.73
## 5 Abilene  2000     5 10590000           10.6 
## 6 Abilene  2000     6 13910000           13.9

Ahora vamos a crear la función. Esta función tomaría 2 argumentos, un marco de datos y un nombre sin procesar de una columna de fecha.create_ymq()df

small_df <- txhousing %>% 
  mutate(date = lubridate::as_date(glue::glue("{year}-{month}-01"))) %>% 
  select(city, date, sales, volume)

create_ymq <- function(df, date_col) {
  stopifnot(inherits(df, "data.frame"))
  stopifnot(class(df %>% dplyr::pull({{date_col}})) == 'Date')
  dplyr::mutate(df,
                Year = lubridate::year({{date_col}}),
                nHalf = lubridate::semester({{date_col}}),
                yHalf = lubridate::semester({{date_col}}, with_year = TRUE),
                dHalf = paste0(lubridate::semester({{date_col}}), "H", format({{date_col}},"%y")),
                nQtr = lubridate::quarter({{date_col}}),
                yQtr = lubridate::quarter({{date_col}}, with_year = TRUE),
                dQtr = paste0(lubridate::quarter({{date_col}}),"Q", format({{date_col}},"%y")),
                Month = lubridate::month({{date_col}}),
                yMonth = as.numeric(format({{date_col}}, "%Y.%m")),
                dMonth = format({{date_col}}, "%b %Y")
                )
}

create_ymq(df = small_df, date_col = date) %>% glimpse()
## Rows: 8,602
## Columns: 14
## $ city   <chr> "Abilene", "Abilene", "Abilene", "Abilene", "Abilene", "Abilene~
## $ date   <date> 2000-01-01, 2000-02-01, 2000-03-01, 2000-04-01, 2000-05-01, 20~
## $ sales  <dbl> 72, 98, 130, 98, 141, 156, 152, 131, 104, 101, 100, 92, 75, 112~
## $ volume <dbl> 5380000, 6505000, 9285000, 9730000, 10590000, 13910000, 1263500~
## $ Year   <dbl> 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 2000, 200~
## $ nHalf  <int> 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 2, 2, 2, ~
## $ yHalf  <dbl> 2000.1, 2000.1, 2000.1, 2000.1, 2000.1, 2000.1, 2000.2, 2000.2,~
## $ dHalf  <chr> "1H00", "1H00", "1H00", "1H00", "1H00", "1H00", "2H00", "2H00",~
## $ nQtr   <int> 1, 1, 1, 2, 2, 2, 3, 3, 3, 4, 4, 4, 1, 1, 1, 2, 2, 2, 3, 3, 3, ~
## $ yQtr   <dbl> 2000.1, 2000.1, 2000.1, 2000.2, 2000.2, 2000.2, 2000.3, 2000.3,~
## $ dQtr   <chr> "1Q00", "1Q00", "1Q00", "2Q00", "2Q00", "2Q00", "3Q00", "3Q00",~
## $ Month  <dbl> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 1, 2, 3, 4, 5, 6, 7, 8, ~
## $ yMonth <dbl> 2000.01, 2000.02, 2000.03, 2000.04, 2000.05, 2000.06, 2000.07, ~
## $ dMonth <chr> "ene. 2000", "feb. 2000", "mar. 2000", "abr. 2000", "may. 2000"~

Ahora que tenemos una función que crea varias columnas relacionadas con la fecha, vamos a crear una función que le permite crear tablas de resumen como ventas anuales por ciudad, volúmenes trimestrales por ciudad, etc.

tx_summary <- function(df, grp_col, sum_col) {
  df %>% 
    group_by(city, {{grp_col}}) %>% 
    summarise("total_{{sum_col}}" := sum({{sum_col}}, na.rm = TRUE), .groups = 'drop')
}
small_df_with_date_cols <- small_df %>% create_ymq(date_col = date)

# Annual Sales per city
small_df_with_date_cols %>% tx_summary(grp_col = Year, sum_col = sales)
## # A tibble: 736 x 3
##    city     Year total_sales
##    <chr>   <dbl>       <dbl>
##  1 Abilene  2000        1375
##  2 Abilene  2001        1431
##  3 Abilene  2002        1516
##  4 Abilene  2003        1632
##  5 Abilene  2004        1830
##  6 Abilene  2005        1977
##  7 Abilene  2006        1997
##  8 Abilene  2007        2003
##  9 Abilene  2008        1651
## 10 Abilene  2009        1634
## # ... with 726 more rows
# Half Yearly volumes per city
small_df_with_date_cols %>% tx_summary(grp_col = yHalf, sum_col = volume)
## # A tibble: 1,472 x 3
##    city    yHalf total_volume
##    <chr>   <dbl>        <dbl>
##  1 Abilene 2000.     55400000
##  2 Abilene 2000.     53175000
##  3 Abilene 2001.     55795000
##  4 Abilene 2001.     58570000
##  5 Abilene 2002.     55305000
##  6 Abilene 2002.     63370000
##  7 Abilene 2003.     58175000
##  8 Abilene 2003.     77500000
##  9 Abilene 2004.     74205000
## 10 Abilene 2004.     85465000
## # ... with 1,462 more rows

https://github.com/strengejacke/ggeffects

# install.packages("ggeffects")

library(ggeffects)
library(splines)

data(efc)
fit <- lm(barthtot ~ c12hour + bs(neg_c_7) * c161sex + e42dep, data = efc)

ggpredict(fit, terms = "c12hour")
## # Predicted values of Total score BARTHEL INDEX
## 
## c12hour | Predicted |         95% CI
## ------------------------------------
##       4 |     67.89 | [65.81, 69.96]
##      12 |     67.07 | [65.10, 69.05]
##      22 |     66.06 | [64.19, 67.94]
##      36 |     64.64 | [62.84, 66.45]
##      49 |     63.32 | [61.51, 65.14]
##      70 |     61.20 | [59.22, 63.17]
##     100 |     58.15 | [55.71, 60.60]
##     168 |     51.26 | [47.27, 55.25]
## 
## Adjusted for:
## * neg_c_7 = 11.83
## * c161sex =  1.76
## *  e42dep =  2.93
library(ggplot2)
mydf <- ggpredict(fit, terms = "c12hour")
ggplot(mydf, aes(x, predicted)) +
  geom_line() +
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = .1)

mydf <- ggpredict(fit, terms = "c12hour")
plot(mydf)

ggpredict(fit, terms = c("neg_c_7", "c161sex", "e42dep"))
## # Predicted values of Total score BARTHEL INDEX
## 
## # c161sex = Male
## #  e42dep = [1] independent
## 
## neg_c_7 | Predicted |          95% CI
## -------------------------------------
##       7 |    102.74 | [95.97, 109.51]
##      12 |    102.27 | [97.10, 107.44]
##      17 |     93.79 | [86.96, 100.63]
##      28 |    164.57 | [95.98, 233.17]
## 
## # c161sex = Female
## #  e42dep = [1] independent
## 
## neg_c_7 | Predicted |           95% CI
## --------------------------------------
##       7 |    109.54 | [105.20, 113.87]
##      12 |     99.81 | [ 95.94, 103.68]
##      17 |     94.90 | [ 90.21,  99.60]
##      28 |     90.26 | [ 71.79, 108.74]
## 
## # c161sex = Male
## #  e42dep = [2] slightly dependent
## 
## neg_c_7 | Predicted |          95% CI
## -------------------------------------
##       7 |     83.73 | [77.32,  90.14]
##      12 |     83.26 | [78.95,  87.58]
##      17 |     74.79 | [68.68,  80.89]
##      28 |    145.57 | [77.00, 214.14]
## 
## # c161sex = Female
## #  e42dep = [2] slightly dependent
## 
## neg_c_7 | Predicted |         95% CI
## ------------------------------------
##       7 |     90.53 | [86.71, 94.35]
##      12 |     80.80 | [78.17, 83.44]
##      17 |     75.90 | [72.29, 79.51]
##      28 |     71.26 | [53.07, 89.45]
## 
## # c161sex = Male
## #  e42dep = [3] moderately dependent
## 
## neg_c_7 | Predicted |          95% CI
## -------------------------------------
##       7 |     64.72 | [58.28,  71.16]
##      12 |     64.26 | [60.30,  68.21]
##      17 |     55.78 | [50.04,  61.52]
##      28 |    126.56 | [57.98, 195.14]
## 
## # c161sex = Female
## #  e42dep = [3] moderately dependent
## 
## neg_c_7 | Predicted |         95% CI
## ------------------------------------
##       7 |     71.52 | [67.59, 75.45]
##      12 |     61.79 | [59.79, 63.80]
##      17 |     56.89 | [53.86, 59.91]
##      28 |     52.25 | [34.21, 70.29]
## 
## # c161sex = Male
## #  e42dep = [4] severely dependent
## 
## neg_c_7 | Predicted |          95% CI
## -------------------------------------
##       7 |     45.72 | [38.86,  52.57]
##      12 |     45.25 | [41.03,  49.47]
##      17 |     36.77 | [30.97,  42.58]
##      28 |    107.55 | [38.93, 176.18]
## 
## # c161sex = Female
## #  e42dep = [4] severely dependent
## 
## neg_c_7 | Predicted |         95% CI
## ------------------------------------
##       7 |     52.51 | [47.88, 57.15]
##      12 |     42.79 | [40.29, 45.28]
##      17 |     37.88 | [34.66, 41.10]
##      28 |     33.24 | [15.21, 51.28]
## 
## Adjusted for:
## * c12hour = 42.10
# install.packages("prettydoc")
library(prettydoc)
mtcars #se aplica el df_print: paged
##                      mpg cyl  disp  hp drat    wt  qsec vs am gear carb
## Mazda RX4           21.0   6 160.0 110 3.90 2.620 16.46  0  1    4    4
## Mazda RX4 Wag       21.0   6 160.0 110 3.90 2.875 17.02  0  1    4    4
## Datsun 710          22.8   4 108.0  93 3.85 2.320 18.61  1  1    4    1
## Hornet 4 Drive      21.4   6 258.0 110 3.08 3.215 19.44  1  0    3    1
## Hornet Sportabout   18.7   8 360.0 175 3.15 3.440 17.02  0  0    3    2
## Valiant             18.1   6 225.0 105 2.76 3.460 20.22  1  0    3    1
## Duster 360          14.3   8 360.0 245 3.21 3.570 15.84  0  0    3    4
## Merc 240D           24.4   4 146.7  62 3.69 3.190 20.00  1  0    4    2
## Merc 230            22.8   4 140.8  95 3.92 3.150 22.90  1  0    4    2
## Merc 280            19.2   6 167.6 123 3.92 3.440 18.30  1  0    4    4
## Merc 280C           17.8   6 167.6 123 3.92 3.440 18.90  1  0    4    4
## Merc 450SE          16.4   8 275.8 180 3.07 4.070 17.40  0  0    3    3
## Merc 450SL          17.3   8 275.8 180 3.07 3.730 17.60  0  0    3    3
## Merc 450SLC         15.2   8 275.8 180 3.07 3.780 18.00  0  0    3    3
## Cadillac Fleetwood  10.4   8 472.0 205 2.93 5.250 17.98  0  0    3    4
## Lincoln Continental 10.4   8 460.0 215 3.00 5.424 17.82  0  0    3    4
## Chrysler Imperial   14.7   8 440.0 230 3.23 5.345 17.42  0  0    3    4
## Fiat 128            32.4   4  78.7  66 4.08 2.200 19.47  1  1    4    1
## Honda Civic         30.4   4  75.7  52 4.93 1.615 18.52  1  1    4    2
## Toyota Corolla      33.9   4  71.1  65 4.22 1.835 19.90  1  1    4    1
## Toyota Corona       21.5   4 120.1  97 3.70 2.465 20.01  1  0    3    1
## Dodge Challenger    15.5   8 318.0 150 2.76 3.520 16.87  0  0    3    2
## AMC Javelin         15.2   8 304.0 150 3.15 3.435 17.30  0  0    3    2
## Camaro Z28          13.3   8 350.0 245 3.73 3.840 15.41  0  0    3    4
## Pontiac Firebird    19.2   8 400.0 175 3.08 3.845 17.05  0  0    3    2
## Fiat X1-9           27.3   4  79.0  66 4.08 1.935 18.90  1  1    4    1
## Porsche 914-2       26.0   4 120.3  91 4.43 2.140 16.70  0  1    5    2
## Lotus Europa        30.4   4  95.1 113 3.77 1.513 16.90  1  1    5    2
## Ford Pantera L      15.8   8 351.0 264 4.22 3.170 14.50  0  1    5    4
## Ferrari Dino        19.7   6 145.0 175 3.62 2.770 15.50  0  1    5    6
## Maserati Bora       15.0   8 301.0 335 3.54 3.570 14.60  0  1    5    8
## Volvo 142E          21.4   4 121.0 109 4.11 2.780 18.60  1  1    4    2

REGRESION LINEAL

https://www.rpubs.com/Joaquin_AR/223351

La regresión lineal simple consiste en generar un modelo de regresión (ecuación de una recta) que permita explicar la relación lineal que existe entre dos variables. A la variable dependiente o respuesta se le identifica como Y y a la variable predictora o independiente como X.

La varianza del error σ2 se estima a partir del Residual Standar Error (RSE), que puede entenderse como la diferencia promedio que se desvía la variable respuesta de la verdadera línea de regresión.

Cuanto mayor es el sumatorio del cuadrado de los residuos menor la precisión con la que el modelo puede predecir el valor de la variable dependiente a partir de la variable predictora. Los residuos son muy importantes puesto que en ellos se basan las diferentes medidas de la bondad de ajuste del modelo

Un analista de deportes quiere saber si existe una relación entre el número de bateos que realiza un equipo de béisbol y el número de runs que consigue. En caso de existir y de establecer un modelo, podría predecir el resultado del partido.

equipos <- c("Texas","Boston","Detroit","Kansas","St.","New_S.","New_Y.",
             "Milwaukee","Colorado","Houston","Baltimore","Los_An.","Chicago",
             "Cincinnati","Los_P.","Philadelphia","Chicago","Cleveland","Arizona",
             "Toronto","Minnesota","Florida","Pittsburgh","Oakland","Tampa",
             "Atlanta","Washington","San.F","San.I","Seattle")
numero_bateos <- c(5659,  5710, 5563, 5672, 5532, 5600, 5518, 5447, 5544, 5598,
                   5585, 5436, 5549, 5612, 5513, 5579, 5502, 5509, 5421, 5559,
                   5487, 5508, 5421, 5452, 5436, 5528, 5441, 5486, 5417, 5421)
runs <- c(855, 875, 787, 730, 762, 718, 867, 721, 735, 615, 708, 644, 654, 735,
          667, 713, 654, 704, 731, 743, 619, 625, 610, 645, 707, 641, 624, 570,
          593, 556)
datos <- data.frame(equipos,numero_bateos,runs)
head(datos)
##   equipos numero_bateos runs
## 1   Texas          5659  855
## 2  Boston          5710  875
## 3 Detroit          5563  787
## 4  Kansas          5672  730
## 5     St.          5532  762
## 6  New_S.          5600  718

cor.test(x = datos$numero_bateos, y = datos$runs, method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  datos$numero_bateos and datos$runs
## t = 4.0801, df = 28, p-value = 0.0003388
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.3209675 0.7958231
## sample estimates:
##      cor 
## 0.610627

El gráfico y el test de correlación muestran una relación lineal, de intensidad considerable (r = 0.61) y significativa (p-value = 0.0003388). Tiene sentido intentar generar un modelo de regresión lineal que permita predecir el número de runs en función del número de bateos del equipo.

0.6.0.2 Cálculo del modelo de regresión lineal simple

modelo_lineal <- lm(runs ~ numero_bateos, datos)
# lm() devuelve el valor de la variable y para x=0 (intersección) junto 
# con la pendiente de la recta.
# Para ver la información del modelo se requiere summary().
summary(modelo_lineal)
## 
## Call:
## lm(formula = runs ~ numero_bateos, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -125.58  -47.05  -16.59   54.40  176.87 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   -2789.2429   853.6957  -3.267 0.002871 ** 
## numero_bateos     0.6305     0.1545   4.080 0.000339 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 66.47 on 28 degrees of freedom
## Multiple R-squared:  0.3729, Adjusted R-squared:  0.3505 
## F-statistic: 16.65 on 1 and 28 DF,  p-value: 0.0003388

La primera columna (Estimate) devuelve el valor estimado para los dos parámetros de la ecuación del modelo lineal (β0 y β1) que equivalen a la ordenada en el origen y la pendiente.

Se muestran los errores estándar, el valor del estadístico t y el p-value (dos colas) de cada uno de los dos parámetros. Esto permite determinar si los parámetros son significativamente distintos de 0, es decir, que tienen importancia en el modelo. En los modelos de regresión lineal simple, el parámetro más informativo suele ser la pendiente. Para el modelo generado, tanto la ordenada en el origen como la pendiente son significativas (p-values < 0.05).

El valor de R2 indica que el modelo calculado explica el 37.29% de la variabilidad presente en la variable respuesta (runs) mediante la variable independiente (número de bateos).

El p-value obtenido en el test F (0.0003388) determina que sí es significativamente superior la varianza explicada por el modelo en comparación a la varianza total. Es el parámetro que determina si el modelo es significativo y por lo tanto se puede aceptar.

El modelo lineal generado sigue la ecuación runs = -2789.2429 + 0.6305 bateos. Por cada unidad que se incrementa el número de bateos, el número de runs aumenta en promedio 0.6305 unidades.

Intervalos de confianza para los parámetros del modelo

confint(modelo_lineal)
##                       2.5 %        97.5 %
## (Intercept)   -4537.9592982 -1040.5264727
## numero_bateos     0.3139863     0.9471137

Representación gráfica del modelo

Además de la línea de mínimos cuadrados es recomendable incluir los límites superior e inferior del intervalo de confianza. Esto permite identificar la región en la que, según el modelo generado y para un determinado nivel de confianza, se encuentra el valor promedio de la variable dependiente.

Para poder representar el intervalo de confianza a lo largo de todo el modelo se recurre a la función predict() para predecir valores que abarquen todo el eje X. Se añaden al gráfico líneas formadas por los límites superiores e inferiores calculados para cada predicción.

# Se genera una secuencia de valores x_i que abarquen todo el rango de las
# observaciones de la variable X
puntos <- seq(from = min(datos$numero_bateos),
              to = max(datos$numero_bateos),
              length.out = 100)
# Se predice el valor de la variable Y junto con su intervalo de confianza para
# cada uno de los puntos generados. En la función predict() hay que nombrar a 
# los nuevos puntos con el mismo nombre que la variable X del modelo.
# Devuelve una matriz.
limites_intervalo <- predict(object = modelo_lineal,
                             newdata = data.frame(numero_bateos = puntos),
                             interval = "confidence", level = 0.95)
head(limites_intervalo, 3)
##        fit      lwr      upr
## 1 626.4464 584.5579 668.3350
## 2 628.3126 587.1743 669.4509
## 3 630.1788 589.7830 670.5745

función geom_smooth() del paquete ggplot2 genera la regresión y su intervalo de forma directa

Por defecto incluye la región de 95% de confianza.

Verificar condiciones para poder aceptar un modelo lineal

Relación lineal entre variable dependiente e independiente:

Se calculan los residuos para cada observación y se representan (scatterplot). Si las observaciones siguen la línea del modelo, los residuos se deben distribuir aleatoriamente entorno al valor 0.

# La función lm() calcula y almacena los valores predichos por el modelo y los residuos.
# 
datos$prediccion <- modelo_lineal$fitted.values
datos$residuos   <- modelo_lineal$residuals
head(datos)
##   equipos numero_bateos runs prediccion  residuos
## 1   Texas          5659  855   779.0395  75.96048
## 2  Boston          5710  875   811.1976  63.80243
## 3 Detroit          5563  787   718.5067  68.49328
## 4  Kansas          5672  730   787.2367 -57.23667
## 5     St.          5532  762   698.9597  63.04033
## 6  New_S.          5600  718   741.8371 -23.83707

Los residuos se distribuyen de forma aleatoria entorno al 0 por lo que se acepta la linealidad.

Distribución normal de los residuos:

Los residuos se deben distribuir de forma normal con media 0. Para comprobarlo se recurre a histogramas, a los cuantiles normales o a un test de contraste de normalidad.

shapiro.test(modelo_lineal$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  modelo_lineal$residuals
## W = 0.96144, p-value = 0.337

Tanto la representación gráfica como el contraste de hipótesis confirman la distribución normal de los residuos.

Varianza constante de los residuos (Homocedasticidad):

La variabilidad de los residuos debe de ser constante a lo largo del eje X. Un patrón cónico es indicativo de falta de homogeneidad en la varianza.

# Test de Breush-Pagan
# install.packages("lmtest")
library(lmtest)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
bptest(modelo_lineal)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo_lineal
## BP = 0.01269, df = 1, p-value = 0.9103

Ni la representación gráfica ni el contraste de hipótesis muestran evidencias que haga sospechar falta de homocedasticidad.

Autocorrelación de residuos:

Cuando se trabaja con intervalos de tiempo, es muy importante comprobar que no existe aoutocorrelación de los residuos, es decir que son independientes. Esto puede hacerse detectando visualmente patrones en la distribución de los residuos cuando se ordenan según se han registrado o con el test de Durbin-Watson dwt() del paquete Car.

Identificación de valores atípicos: outliers, leverage y observaciones influyentes

Outlier u observación atípica: Observaciones que no se ajustan bien al modelo. El valor real se aleja mucho del valor predicho, por lo que su residuo es excesivamente grande. En una representación bidimensional se corresponde con desviaciones en el eje Y.

Observación influyente: Observación que influye sustancialmente en el modelo, su exclusión afecta al ajuste. No todos los outliers tienen por qué ser influyentes.

Observación con alto leverage: Observación con un valor extremo para alguno de los predictores. En una representación bidimensional se corresponde con desviaciones en el eje X. Son potencialmente puntos influyentes.

Independientemente de que el modelo se haya podido aceptar, siempre es conveniente identificar si hay algún posible outlier, observación con alto leverage u observación altamente influyente, puesto que podría estar condicionando en gran medida el modelo. La eliminación de este tipo de observaciones debe de analizarse con detalle y dependiendo de la finalidad del modelo. Si el fin es predictivo, un modelo sin estas observaciones puede lograr mayor precisión la mayoría de casos. Sin embargo, es muy importante prestar atención a estos valores ya que, de no ser errores de medida, pueden ser los casos más interesantes. El modo adecuado a proceder cuando se sospecha de algún posible valor atípico o influyente es calcular el modelo de regresión incluyendo y excluyendo dicho valor.

datos %>% filter(abs(studentized_residual) > 3)
##   equipos numero_bateos runs prediccion residuos studentized_residual
## 1  New_Y.          5518  867    690.132  176.868             3.092876
which(abs(datos$studentized_residual) > 3)
## [1] 7

El estudio de los residuos studentized identifica al equipo de New_Y. como una posible observación atípica. Esta observación ocupa la posición 7 en la tabla de datos.

El hecho de que un valor sea atípico o con alto grado de leverage no implica que sea influyente en el conjunto del modelo. Sin embargo, si un valor es influyente, suele ser o atípico o de alto leverage. En R se dispone de la función outlierTest() del paquete car y de las funciones influence.measures(), influencePlot() y hatvalues() para identificar las observaciones más influyentes en el modelo.

library(car)
## Loading required package: carData
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
## The following object is masked from 'package:purrr':
## 
##     some
summary(influence.measures(model = modelo_lineal))
## Potentially influential observations of
##   lm(formula = runs ~ numero_bateos, data = datos) :
## 
##   dfb.1_ dfb.nmr_ dffit cov.r   cook.d hat    
## 2 -0.53   0.54     0.58  1.27_*  0.17   0.22_*
## 7  0.05  -0.04     0.58  0.61_*  0.13   0.03

      StudRes        Hat      CookD
2   1.0914283 0.22133381 0.16815163
4  -0.9331751 0.15252728 0.07872749
7   3.0928757 0.03349684 0.12693385
10 -2.0622189 0.06333282 0.12881098

Las funciones influence.measures() e influencePlot() detectan la observación 7 como atípica pero no significativamente influyente. Sí detectan como influyente la observación que ocupa la segunda posición. Para evaluar hasta qué punto condiciona el modelo, se recalcula la recta de mínimos cuadrados excluyendo esta observación.

La eliminación del valor identificado como influyente apenas cambia la recta de mínimos cuadrados. Para conocer con exactitud el resultado de excluir la observación se comparan las pendientes de ambos modelos.

lm(formula = runs ~ numero_bateos, data = datos)$coefficients
##   (Intercept) numero_bateos 
##   -2789.24289       0.63055
lm(formula = runs ~ numero_bateos, data = datos[-2,])$coefficients
##   (Intercept) numero_bateos 
## -2335.7478247     0.5479527

Conclusión

Dado que se satisfacen todas las condiciones para considerar válido un modelo de regresión lineal por mínimos cuadrados y que el p-value indica que el ajuste es significativo, se puede aceptar el modelo lineal. A pesar de ello, el valor de R2 no es muy alto por lo que el número de bateos no es muy buen predictor del número de runs.

Bibliografía

Linear Models with R by Julian J.Faraway libro An Introduction to Statistical Learning: with Applications in R (Springer Texts in Statistics) libro OpenIntro Statistics: Fourth Edition by David Diez, Mine Çetinkaya-Rundel, Christopher Barr libro Extending the Linear Model with R: Generalized Linear, Mixed Effects and Nonparametric Regression Models by Julian J.Faraway libro Introduction to Machine Learning with Python: A Guide for Data Scientists libro Points of Significance: Association, correlation and causation. Naomi Altman & Martin Krzywinski Nature Methods Points of Significance: Simple linear regression Naomi Altman & Martin Krzywinski. Nature Methods Resampling Data: Using a Statistical Jackknife S. Sawyer | Washington University | March 11, 2005 https://en.wikipedia.org/wiki/Resampling_(statistics)#Jackknife The Trusty Jackknife Method identifies outliers and bias in statistical estimates by I. Elaine Allen and Christopher A. Seaman

sesion_info <- devtools::session_info()
dplyr::select(
  tibble::as_tibble(sesion_info$packages),
  c(package, loadedversion, source)
)
## # A tibble: 161 x 3
##    package    loadedversion source        
##    <chr>      <chr>         <chr>         
##  1 abind      1.4-5         CRAN (R 4.1.1)
##  2 assertthat 0.2.1         CRAN (R 4.1.2)
##  3 backports  1.4.1         CRAN (R 4.1.2)
##  4 base64enc  0.1-3         CRAN (R 4.1.1)
##  5 bitops     1.0-7         CRAN (R 4.1.1)
##  6 broom      0.7.11        CRAN (R 4.1.2)
##  7 bslib      0.3.1         CRAN (R 4.1.2)
##  8 cachem     1.0.6         CRAN (R 4.1.2)
##  9 callr      3.7.0         CRAN (R 4.1.2)
## 10 car        3.0-12        CRAN (R 4.1.2)
## # ... with 151 more rows

LEER ARCHIVOS XML

https://analisisydecision.es/trucos-r-leer-archivos-xml-con-r/

# Install packages
# install.packages("XML")
# install.packages("curl")

# Load packages
library(XML)
## Warning: package 'XML' was built under R version 3.6.3
library(curl)

  1. This is a limerick written by Claus Ekstrøm: https://yihui.org/en/2018/06/xaringan-math-limerick/.↩︎